Modelo Basado en Agentes

Author

Domingo Martínez

Code
# Lista de paquetes necesarios
required_packages <- c("tidyverse",
                       "NetLogoR",
                       "SpaDES.core")

# Función para instalar y cargar paquetes automáticamente
load_or_install <- function(packages) {
  for (pkg in packages) {
    if (!require(pkg, character.only = TRUE)) {
      install.packages(pkg, dependencies = TRUE)
      library(pkg, character.only = TRUE)
    }
  }
}

# Llamada a la función con los paquetes necesarios
load_or_install(required_packages)
Cargando paquete requerido: tidyverse
── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.4     ✔ readr     2.1.5
✔ forcats   1.0.0     ✔ stringr   1.5.1
✔ ggplot2   3.5.1     ✔ tibble    3.2.1
✔ lubridate 1.9.4     ✔ tidyr     1.3.1
✔ purrr     1.0.2     
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
Cargando paquete requerido: NetLogoR

Cargando paquete requerido: SpaDES.core

Cargando paquete requerido: quickPlot

Cargando paquete requerido: reproducible


Adjuntando el paquete: 'SpaDES.core'


The following object is masked from 'package:NetLogoR':

    updateList


The following objects are masked from 'package:lubridate':

    dmonths, dweeks, dyears


The following objects are masked from 'package:stats':

    end, start


The following object is masked from 'package:utils':

    citation

1 Ejemplo COVID-19

Code
nTime1=20 # Vamos a simular 20 díás de intereacción
Infected <- data.frame(Day = 1:nTime1) # un data frame para contar infectados
Susceptible <- data.frame(Day = 1:nTime1) # Un data frame para contar susceptibles
Exposed <- data.frame(Day = 1:nTime1) # Un  data frame para contar expuestos

# Declaramos un ciclo 

for (m in 1:20) {
  #####################################################################
  usan_cubrebocas <- data.frame(AgentNo = 1,
                                    State = "S", 
                                    ObeySocialNorm = 1,
                                    Risk_seeking = rpois(1,2),
                                    stringsAsFactors = F,
                                    Social_profile = "Usa_cubrebocas")
    
  no_usan_cubrebocas <- data.frame(AgentNo = 1,
                                    State = "E",
                                    ObeySocialNorm = 0,
                                    Risk_seeking = rpois(1,2),
                                    stringsAsFactors = F,
                                    Social_profile = "No_usa_cubrebocas")
  #####################################################################
  # Now, lets create a population for each kind of agent
  
  # define size of each population
  Usa_population = 10 # Definimos población que sí usa cubrebocas
  No_usa_population = 1 # Definimos población que NO usa cubrebocas
  length.population = Usa_population + No_usa_population
  # create population Usa cubrebocas
    for (i in 2:Usa_population) {
      usa_cubrebocas <- data.frame(AgentNo = 1,
                                    State = "S", 
                                    ObeySocialNorm = 1,
                                    Risk_seeking = rpois(1,2),
                                    stringsAsFactors = F,
                                    Social_profile = "Usa_cubrebocas")
     
      usan_cubrebocas <- rbind(usa_cubrebocas,usan_cubrebocas)
    }
    
    # create population NO usa cubrebocas
    
    for (i in 2:No_usa_population) {
      no_usa_cubrebocas <- data.frame(AgentNo = 1,
                                    State = "E",
                                    ObeySocialNorm = 0,
                                    Risk_seeking = rpois(1,2),
                                    stringsAsFactors = F,
                                    Social_profile = "No_usa_cubrebocas")
      
      no_usan_cubrebocas <- rbind(no_usan_cubrebocas, no_usa_cubrebocas)
    }
  #################################################################
  # Now we joint three data sets in only one.
  population <- rbind(usan_cubrebocas, no_usa_cubrebocas)
  population <- population[,-1] # Eliminamos la última columna
  glimpse(population) # Revisamos la estructura de la población
  #################################################################
    nTime1 <- 20 # Set how far the model will run.
    # Now we will create a data frame filled with zeros.
    Out1 <- data.frame(Day = 1:nTime1,
                       S = rep(0, nTime1),
                       E = rep(0, nTime1),
                       I = rep(0, nTime1)
                       )
    # Move the people through time.
    for (k in 1:nTime1) {
      # A cycle for moving people through time.
      for (i in 1:length.population) {
        # Here we determine the Risk seeking behavior, i.e. socialize. 
        # In other words, if they like to meet others during the lock down.
        # Number [i] indicates no wearing mask agent number i.
        # Now, we depict how many persons will meet person [i]
        Meet1 <- population$Risk_seeking[ i ] # 
        Meet1 
        # Now, let's grab the agents they will meet
        Meet2 <- sample(1:length.population, Meet1, replace = T)
        # Add the subject i to the group
        Meet2 <- c(Meet2,i)
        Meet2
        # Check if someone in the group has been exposed.
        exposed.table <- table(population[Meet2,1] == "E")
        true_values <- exposed.table["TRUE"]
        true_values[is.na(true_values)] <- 0
        true_values <- as.numeric(true_values)
        # If someone got exposed, then, only "S" people in group change status to "E"
        if (true_values > 0) {
          # Grab only agents who are susceptible from Meet2
          a <- population[Meet2,1]
          a <-(a == "S")
          b <- data.frame(A = Meet2,
                          B = a)
          b
          Meet2.1 <- subset(b,B == TRUE, select = A)
          Meet2.1
          population$State[Meet2.1$A] <- "E"
          
         # Check if all people into the group obey social norm 
          obey.table <- table(population[Meet2,2] == 0)
          true_values2 <- obey.table["TRUE"]
          true_values2[is.na(true_values2)] <- 0
          true_values2 <- as.numeric(true_values2)
          # If, at least, one person don't obey the social norm, then the exposed people can get infected, under certain probability. 
          
            if (true_values2 > 0) {
            # Grab people into Meet2
              for (j in Meet2) {
              # Let's assign a probability to get infected to each subject into Meet2.
                prob.to.get.infected <- runif(1,0,1)
                if (prob.to.get.infected > 0.7) {
                  population$State[j] <- "I"
                }
              }
            }
         }
      }
      # Now we record the results progressively.
      Out1$S[k]<-length(population$State[population$State == "S"])
      Out1$E[k]<-length(population$State[population$State == "E"])
      Out1$I[k]<-length(population$State[population$State == "I"])
      
    }
    
    Out1
  #################################################################
  Infected<-cbind(Infected,Out1$I)
  Susceptible<-cbind(Susceptible,Out1$S)
  Exposed<-cbind(Exposed,Out1$E)
  #################################################################
}
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 1, 0, 3, 2, 1, 0, 2, 1, 2, 1, 2
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 2, 4, 4, 1, 2, 3, 3, 1, 2, 3, 1
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 0, 2, 3, 3, 3, 1, 1, 4, 1, 0, 3
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 3, 3, 3, 1, 1, 5, 3, 1, 3, 3, 2
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 1, 1, 3, 1, 1, 3, 1, 3, 0, 1, 1
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 1, 0, 2, 2, 1, 3, 1, 1, 0, 1, 3
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 2, 1, 3, 4, 0, 1, 2, 1, 3, 6, 3
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 0, 2, 1, 1, 5, 2, 3, 3, 1, 0, 3
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 3, 3, 2, 1, 2, 2, 4, 1, 3, 0, 2
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 2, 3, 0, 4, 3, 4, 1, 2, 2, 0, 2
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 4, 4, 1, 0, 0, 1, 3, 1, 0, 2, 5
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 3, 4, 1, 4, 1, 2, 3, 1, 3, 2, 1
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 2, 1, 0, 3, 2, 4, 7, 2, 2, 3, 1
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 4, 4, 3, 2, 3, 1, 2, 1, 1, 1, 0
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 3, 2, 4, 3, 2, 4, 2, 2, 3, 1, 1
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 1, 3, 1, 1, 3, 2, 1, 1, 0, 0, 1
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 4, 3, 0, 3, 1, 5, 1, 1, 2, 3, 2
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 2, 1, 1, 1, 1, 0, 2, 1, 3, 0, 0
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 2, 1, 1, 3, 2, 6, 1, 3, 2, 2, 6
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Rows: 11
Columns: 4
$ State          <chr> "S", "S", "S", "S", "S", "S", "S", "S", "S", "S", "E"
$ ObeySocialNorm <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0
$ Risk_seeking   <int> 3, 2, 0, 4, 2, 2, 2, 3, 2, 0, 2
$ Social_profile <chr> "Usa_cubrebocas", "Usa_cubrebocas", "Usa_cubrebocas", "…
Code
# Now, let’s get the average of infected people every day.
Infected<-apply(Infected[,-1],1, mean)
Susceptible<-apply(Susceptible[,-1],1, mean)
Exposed<-apply(Exposed[,-1],1, mean)
nTime1=20
final.results <- data.frame(Day = 1:nTime1,
                            Infected,
                            Susceptible,
                            Exposed)

2 Gráfico Dinámico COVID-19

Code
# Let's take a look
library(plotly)

Adjuntando el paquete: 'plotly'
The following object is masked from 'package:ggplot2':

    last_plot
The following object is masked from 'package:stats':

    filter
The following object is masked from 'package:graphics':

    layout
Code
fig1 <- plot_ly(final.results)%>%
  
  add_trace(type = 'scatter', mode = 'lines', x = ~Day, y = ~Exposed, hovertemplate = 'Exposed: %{y:}<extra></extra>', name = "Exposed",line = list(color = 'rgb(105,105,105)', width = 2, dash = 'dash'))%>%
  
  add_trace(type = 'scatter', mode = 'lines', x = ~Day, y = ~Infected, hovertemplate = 'Infected: %{y:}<extra></extra>', name = "Infected", line = list(color = 'rgb(105,105,105)', width = 2))%>%
  
  add_trace(type = 'scatter', mode = 'lines',x = ~Day, y = ~Susceptible, hovertemplate = 'Susceptible: %{y:}<extra></extra>', name = "Susceptible", line = list(color = 'rgb(128, 128, 128)', width = 2, dash = 'dot') )%>%
  layout(showlegend = T)
  

fig1 <- fig1 %>%

  layout(hovermode="x unified",

         xaxis = list(zerolinecolor = '#ffff',

                      zerolinewidth = 2,

                      gridcolor = 'ffff', 
                      
                      title = "Day"),

         yaxis = list(zerolinecolor = '#ffff',

                      zerolinewidth = 2,

                      gridcolor = 'ffff', 
                      
                      title = "Number of people"),
         
         title = 'COVID-19 propagation',

         plot_bgcolor='#DCDCDC')
fig1
Code
library(htmlwidgets)
saveWidget(fig1, "COVID-19 propagation.html", selfcontained = F, libdir = "lib")

3 Gráfico Estático COVID-19

Code
ggplot(data = final.results, mapping = aes(x = Day)) +
  geom_line(aes(y = Infected, color = "Infected", linetype = "Infected")) +
  geom_line(aes(y = Exposed, color = "Exposed", linetype = "Exposed")) +
  geom_line(aes(y = Susceptible, color = "Susceptible", linetype = "Susceptible")) +
  scale_color_manual(
    values = c("Infected" = "black", "Exposed" = "gray40", "Susceptible" = "gray60")
  ) +
  scale_linetype_manual(
    values = c("Infected" = "solid", "Exposed" = "dashed", "Susceptible" = "dotted")
  ) +
  labs(
    title = "COVID-19 propagation",
    x = "Day",
    y = "Number of people",
    color = "State",
    linetype = "Satate"
  ) +
  theme_minimal() +
  guides(color = guide_legend(title = "State"), linetype = guide_legend(title = "State"))

Code
ggsave("Fig_01_COVID-19_propagation.png", dpi = "retina", width = 12, height = 8, units = "cm")

4 Ejemplo depredador-presa con NetLogoR

Code
library(NetLogoR)

# Crear un mundo de 50x50 celdas
world <- createWorld(minPxcor = 0, maxPxcor = 50, minPycor = 0, maxPycor = 50)

# Crear 50 conejos y 10 lobos
conejos <- createTurtles(n = 50, world = world)
lobos <- createTurtles(n = 10, world = world)

# Simulación de movimiento aleatorio
move_agents <- function(agents) {
  agents <- fd(agents, dist = 1)  # Mover 1 unidad
  agents <- right(agents, angle = runif(nrow(agents), -45, 45))  # Giro aleatorio
  return(agents)
}

# Ejecutar 10 pasos de simulación
for (i in 1:10) {
  conejos <- move_agents(conejos)
  lobos <- move_agents(lobos)
}
print(conejos)
   who    heading    prevX    prevY  breed   color
1    0   8.198790 26.24847 33.12423 turtle #FF0000
2    1  43.344252 33.25555 25.37099 turtle #FF1F00
3    2 162.725163 19.96019 25.00817 turtle #FF3D00
4    3 124.261710 28.88568 19.91214 turtle #FF5C00
5    4  41.853179 22.65584 32.55115 turtle #FF7A00
6    5 262.197367 17.41328 22.25541 turtle #FF9900
7    6  90.549191 26.60008 17.03939 turtle #FFB800
8    7 268.945258 23.05904 17.61324 turtle #FFD600
9    8  86.535210 31.51669 19.40970 turtle #FFF500
10   9 277.990557 22.27133 18.96763 turtle #EBFF00
11  10   8.593784 30.00088 31.10763 turtle #CCFF00
12  11 296.302370 23.76303 19.79960 turtle #ADFF00
13  12 266.984047 19.68802 19.89843 turtle #8FFF00
14  13 150.146754 32.03693 20.91713 turtle #70FF00
15  14  41.635784 31.41335 30.52465 turtle #52FF00
16  15 113.092109 23.74508 31.77433 turtle #33FF00
17  16 192.303674 27.62262 17.46468 turtle #14FF00
18  17 325.879074 23.03507 32.39918 turtle #00FF0A
19  18  53.986845 30.28605 31.57751 turtle #00FF29
20  19 259.753447 16.84464 25.53777 turtle #00FF47
21  20  37.793022 32.70817 28.23656 turtle #00FF66
22  21 132.452182 33.21535 26.43084 turtle #00FF85
23  22  18.913920 21.31006 32.84683 turtle #00FFA3
24  23 321.006541 20.32500 31.01116 turtle #00FFC2
25  24 319.008124 22.70303 30.19966 turtle #00FFE0
26  25  14.438925 26.15073 33.39778 turtle #00FFFF
27  26 299.075161 17.04795 24.72501 turtle #00E0FF
28  27 284.247656 16.53820 25.98567 turtle #00C2FF
29  28  71.211736 31.23993 30.32171 turtle #00A3FF
30  29  18.619206 25.93608 33.11190 turtle #0085FF
31  30  56.335889 32.90751 23.14173 turtle #0066FF
32  31 226.845759 17.66557 28.13881 turtle #0047FF
33  32 336.393683 31.15382 30.62223 turtle #0029FF
34  33 211.833371 25.99520 17.99879 turtle #000AFF
35  34  15.329774 24.78026 33.15767 turtle #1400FF
36  35 133.683590 27.20574 16.96789 turtle #3300FF
37  36  23.241553 31.41403 30.88662 turtle #5200FF
38  37 154.231130 30.28627 19.40791 turtle #7000FF
39  38 158.296627 22.27627 16.77512 turtle #8F00FF
40  39  50.411272 30.09357 29.15827 turtle #AD00FF
41  40  77.019637 32.10273 30.16433 turtle #CC00FF
42  41 194.184365 18.77066 21.40319 turtle #EB00FF
43  42 172.187758 19.42394 22.16570 turtle #FF00F5
44  43 259.702062 18.78559 26.65913 turtle #FF00D6
45  44 173.517672 30.06259 19.66182 turtle #FF00B8
46  45  91.228362 25.32259 19.11826 turtle #FF0099
47  46 167.242317 33.11843 22.29189 turtle #FF007A
48  47 206.885950 16.78295 25.19494 turtle #FF005C
49  48  62.801266 29.38520 32.02261 turtle #FF003D
50  49 131.045181 26.17241 32.10428 turtle #FF001F
Code
print(lobos)
   who   heading    prevX    prevY  breed   color
1    0  12.69352 27.71302 29.59000 turtle #FF0000
2    1 190.37016 30.83267 27.16834 turtle #FF9900
3    2 300.27424 17.62165 29.22931 turtle #CCFF00
4    3  10.77593 29.46681 31.25102 turtle #33FF00
5    4 310.62729 18.13995 25.78771 turtle #00FF66
6    5 296.02763 17.32110 25.55525 turtle #00FFFF
7    6 211.59396 18.48857 23.74900 turtle #0066FF
8    7  92.79403 31.91460 28.22419 turtle #3300FF
9    8  94.33338 26.77143 18.56701 turtle #CC00FF
10   9 331.67905 16.89678 26.89450 turtle #FF0099

5 Ejemplo Dinámica de población con SpaDES

Code
library(SpaDES.core)

# Definir un módulo de simulación simple
params <- list(
  growthRate = 0.1,
  deathRate = 0.05
)

simulate_population <- function(steps, init_pop = 100) {
  pop <- init_pop
  for (t in 1:steps) {
    births <- rbinom(1, pop, params$growthRate)
    deaths <- rbinom(1, pop, params$deathRate)
    pop <- pop + births - deaths
    cat("Paso", t, "- Población:", pop, "\n")
  }
}

simulate_population(steps = 10)
Paso 1 - Población: 104 
Paso 2 - Población: 110 
Paso 3 - Población: 109 
Paso 4 - Población: 109 
Paso 5 - Población: 114 
Paso 6 - Población: 117 
Paso 7 - Población: 124 
Paso 8 - Población: 132 
Paso 9 - Población: 138 
Paso 10 - Población: 151